{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "# utility modules\n",
    "import glob\n",
    "import os\n",
    "import sys\n",
    "import re\n",
    "import io\n",
    "\n",
    "# the usual suspects\n",
    "import numpy as np\n",
    "import cartopy.crs as ccrs\n",
    "from cartopy.io import shapereader\n",
    "from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER\n",
    "import cartopy.io.img_tiles as cimgt\n",
    "import matplotlib.pyplot as plt\n",
    "import matplotlib.cm as cm\n",
    "# import matplotlib.image as mpimg\n",
    "\n",
    "# modules you'll need if you're downloading the data\n",
    "from icepyx import icesat2data as ipd\n",
    "import shutil\n",
    "import geopandas as gpd\n",
    "import pandas as pd\n",
    "from shapely.geometry import Point, Polygon\n",
    "\n",
    "import fiona\n",
    "import pyproj\n",
    "import h5py\n",
    "\n",
    "# # run matplotlib in 'widget' mode\n",
    "# %matplotlib widget\n",
    "# %load_ext autoreload\n",
    "# %autoreload 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "from readers.read_HDF5_ATL13 import read_HDF5_ATL13\n",
    "ATL13_file = 'ATL13_20190128160153_04770201_003_01.h5'\n",
    "IS2_atl13_mds, IS2_atl13_attrs, IS2_atl13_beams = read_HDF5_ATL13(ATL13_file)\n",
    "beam = ['gt1l','gt1r','gt2l','gt2r','gt3l','gt3r']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "4422.655192057292\n",
      "0.37747750118800943\n"
     ]
    }
   ],
   "source": [
    "fiona.drvsupport.supported_drivers['kml'] = 'rw' \n",
    "fiona.drvsupport.supported_drivers['KML'] = 'rw'\n",
    "lake_filepath = os.getcwd()+'/sikkimlakes.kml'\n",
    "gdfl = gpd.read_file(lake_filepath) \n",
    "mask = gdfl.geometry.unary_union\n",
    " \n",
    "D13 = {}\n",
    "stdev = []\n",
    "surface = 0\n",
    "for i in beam:\n",
    "    D13[i] = {}    \n",
    "    flag = np.zeros(IS2_atl13_mds[i]['segment_lon'].size, dtype=bool)\n",
    "    for j in range(0, IS2_atl13_mds[i]['segment_lon'].size):\n",
    "        flag[j] = (Point(IS2_atl13_mds[i]['segment_lon'][j],IS2_atl13_mds[i]['segment_lat'][j]).within(mask)) & (IS2_atl13_mds[i]['segment_lat'][j]>28.3)\n",
    "    D13[i]['lon'] = IS2_atl13_mds[i]['segment_lon'][flag]\n",
    "    D13[i]['lat'] = IS2_atl13_mds[i]['segment_lat'][flag]\n",
    "    D13[i]['h_ortho'] = IS2_atl13_mds[i]['ht_ortho'][flag]\n",
    "    D13[i]['ice_flag'] = IS2_atl13_mds[i]['ice_flag'][flag]\n",
    "    stdev = np.append(stdev, D13[i]['h_ortho'][~np.isnan(D13[i]['h_ortho'])])\n",
    "    # calculate the mean surface height of lake\n",
    "    surface = surface + np.average(D13[i]['h_ortho'][~np.isnan(D13[i]['h_ortho'])])\n",
    "surface /= len(beam)\n",
    "print(surface)\n",
    "stdev = np.std(stdev)\n",
    "print(stdev)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX0AAAD4CAYAAAAAczaOAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAASQElEQVR4nO3de4xc51nH8e8Pu3VDiXshTgleBxvJlZpUVdQMViDQkgKNSUsSrjWipBIVFiGgtKQUrEpF/MO9UAJKJKutgiltGi4hUdW0CZTbH0nTcZs0t4YuJDSuU+wKCQwVbp0+/DHH6bAe7856d2fGfr8fabRnnnPe2fdJvL9598zZmVQVkqQ2fMO0JyBJmhxDX5IaYuhLUkMMfUlqiKEvSQ1ZP+0JLOWcc86prVu3TnsaknRa2b9//5eqatPC+syH/tatW+n3+9OehiSdVpL826i6p3ckqSGGviQ1xNCXpIYY+pLUEENfkhoy81fvnIp3veENwP+cUN/wgiv4hb0/P/kJSdKMOONW+icLfICj//kR/nj3TZOdkCTNkDMu9E8W+McdPfJPE5qHJM2eMzD0l/C1I9OegSRNTXuh/w1nT3sGkjQ1Z2DoP3/RvRvO/p4JzUOSZs8ZF/o3fOhDnCz4vXpHUuvOyEs2B8EvSVrojFvpS5JObqzQT/JkkoeSPJCkv2Df25JUknOGanuSzCd5PMnlQ/WLu8eZT3JjkqxeK5KkpSxnpX9ZVV1UVb3jhSRbgB8APj9UuwDYBVwI7ARuSrKu230zsBvY3t12rmz6kqTlWOnpnT8A3g7UUO0q4NaqOlpVTwDzwI4k5wEbq+reqipgH3D1Cr+/JGkZxg39Au5Osj/JboAkVwJfqKoHFxy7GXhq6P6Brra5215YP0GS3Un6SfqHDx8ec4qSpKWMe/XOpVV1MMm5wD1JPgu8A3jtiGNHnaevReonFqv2AnsBer3eyGMkScs31kq/qg52Xw8BtwOvBrYBDyZ5EpgDPpXkWxis4LcMDZ8DDnb1uRF1SdKELBn6SZ6f5Ozj2wxW95+sqnOramtVbWUQ6K+sqi8CdwK7kmxIso3BC7b3V9XTwJEkl3RX7VwD3LE2bUmSRhnn9M5LgNu7qyvXAx+oqo+e7OCqeiTJbcCjwDHguqp6ptt9LXALcBZwV3eTJE1IBhfSzK5er1f9fn/pAyVJz0qyf/gS++P8i1xJaoihL0kNMfQlqSGGviQ1xNCXpIYY+pLUEENfkhpi6EtSQwx9SWqIoS9JDTH0Jakhhr4kNcTQl6SGGPqS1BBDX5IaYuhLUkMMfUlqiKEvSQ0x9CWpIYa+JDXE0Jekhhj6ktQQQ1+SGmLoS1JDDH1JaoihL0kNMfQlqSGGviQ1xNCXpIaMFfpJnkzyUJIHkvS72u8m+WySzyS5PckLh47fk2Q+yeNJLh+qX9w9znySG5Nk9VuSJJ3Mclb6l1XVRVXV6+7fA7y8ql4B/DOwByDJBcAu4EJgJ3BTknXdmJuB3cD27rZz5S1IksZ1yqd3quruqjrW3b0PmOu2rwJuraqjVfUEMA/sSHIesLGq7q2qAvYBV69g7pKkZRo39Au4O8n+JLtH7P8Z4K5uezPw1NC+A11tc7e9sH6CJLuT9JP0Dx8+POYUJUlLGTf0L62qVwI/CFyX5FXHdyR5B3AM+LPjpRHja5H6icWqvVXVq6repk2bxpyiJGkpY4V+VR3svh4Cbgd2ACR5E/B64Ke6UzYwWMFvGRo+Bxzs6nMj6pKkCVky9JM8P8nZx7eB1wIPJ9kJ/ApwZVV9eWjIncCuJBuSbGPwgu39VfU0cCTJJd1VO9cAd6xyP5KkRawf45iXALd3V1euBz5QVR9NMg9sAO7p9t1XVT9XVY8kuQ14lMFpn+uq6pnusa4FbgHOYvAawF1IkiYmXz8rM5t6vV71+/1pT0OSTitJ9g9dYv8s/yJXkhpi6EtSQwx9SWqIoS9JDTH0Jakhhr4kNcTQl6SGGPqS1BBDX5IaYuhLUkMMfUlqiKEvSQ0x9CWpIYa+JDXE0Jekhhj6ktQQQ1+SGmLoS1JDDH1JaoihL0kNMfQlqSGGviQ1xNCXpIYY+pLUEENfkhpi6EtSQwx9SWqIoS9JDTH0Jakhhr4kNWSs0E/yZJKHkjyQpN/VXpzkniSf676+aOj4PUnmkzye5PKh+sXd48wnuTFJVr8lSdLJLGelf1lVXVRVve7+rwJ/W1Xbgb/t7pPkAmAXcCGwE7gpybpuzM3AbmB7d9u58hYkSeNayemdq4A/6bb/BLh6qH5rVR2tqieAeWBHkvOAjVV1b1UVsG9ojCRpAsYN/QLuTrI/ye6u9pKqehqg+3puV98MPDU09kBX29xtL6yfIMnuJP0k/cOHD485RUnSUtaPedylVXUwybnAPUk+u8ixo87T1yL1E4tVe4G9AL1eb+QxkqTlG2ulX1UHu6+HgNuBHcC/d6ds6L4e6g4/AGwZGj4HHOzqcyPqkqQJWTL0kzw/ydnHt4HXAg8DdwJv6g57E3BHt30nsCvJhiTbGLxge393CuhIkku6q3auGRojSZqAcU7vvAS4vbu6cj3wgar6aJJPArcleTPweeDHAarqkSS3AY8Cx4DrquqZ7rGuBW4BzgLu6m6SpAnJ4EKa2dXr9arf7097GpJ0Wkmyf+gS+2f5F7mS1BBDX5IaYuhLUkMMfUlqiKEvSQ0x9CWpIYa+JDXE0Jekhhj6ktQQQ1+SGmLoS1JDDH1JaoihL0kNMfQlqSGGviQ1xNCXpIYY+pLUEENfkhpi6EtSQwx9SWqIoS9JDTH0Jakhhr4kNcTQl6SGGPqS1BBDX5IaYuhLUkMMfUlqiKEvSQ0x9CWpIWOHfpJ1ST6d5MPd/YuS3JfkgST9JDuGjt2TZD7J40kuH6pfnOShbt+NSbK67UiSFrOclf71wGND938H+PWqugh4Z3efJBcAu4ALgZ3ATUnWdWNuBnYD27vbzhXNXpK0LGOFfpI54HXAe4bKBWzstl8AHOy2rwJuraqjVfUEMA/sSHIesLGq7q2qAvYBV69CD5KkMa0f87h3A28Hzh6qvQX4WJLfY/Dk8V1dfTNw39BxB7raV7vthfUTJNnN4DcCzj///DGnKElaypIr/SSvBw5V1f4Fu64F3lpVW4C3Au89PmTEw9Qi9ROLVXurqldVvU2bNi01RUnSmMZZ6V8KXJnkCuB5wMYk7wd+iMF5foA/5+unfg4AW4bGzzE49XOg215YlyRNyJIr/araU1VzVbWVwQu0H6+qNzII7Fd3h70G+Fy3fSewK8mGJNsYvGB7f1U9DRxJckl31c41wB2r244kaTHjntMf5WeBP0yyHvhfunPwVfVIktuAR4FjwHVV9Uw35lrgFuAs4K7uJkmakAwupJldvV6v+v3+tKchSaeVJPurqrew7l/kSlJDDH1JaoihL0kNMfQlqSGGviQ1xNCXpIYY+pLUEENfkhpi6EtSQwx9SWqIoS9JDTH0Jakhhr4kNcTQl6SGGPqS1BBDX5IaYuhLUkMMfUlqiKEvSQ0x9CWpIYa+JDXE0Jekhhj6ktQQQ1+SGmLoS1JDDH1JaoihL0kNMfQlqSGGviQ1ZOzQT7IuyaeTfHio9otJHk/ySJLfGarvSTLf7bt8qH5xkoe6fTcmyeq1IklayvplHHs98BiwESDJZcBVwCuq6miSc7v6BcAu4ELgW4G/SfLSqnoGuBnYDdwHfATYCdy1Sr1IkpYw1ko/yRzwOuA9Q+Vrgd+qqqMAVXWoq18F3FpVR6vqCWAe2JHkPGBjVd1bVQXsA65epT4kSWMY9/TOu4G3A18bqr0U+J4kn0jyD0m+o6tvBp4aOu5AV9vcbS+sS5ImZMnQT/J64FBV7V+waz3wIuAS4JeB27pz9KPO09ci9VHfc3eSfpL+4cOHl5qiJGlM45zTvxS4MskVwPOAjUnez2Cl/lfdqZr7k3wNOKerbxkaPwcc7OpzI+onqKq9wF6AXq838olBkrR8S670q2pPVc1V1VYGL9B+vKreCPw18BqAJC8Fngt8CbgT2JVkQ5JtwHbg/qp6GjiS5JLuN4JrgDvWoilJ0mjLuXpnofcB70vyMPAV4E3dqv+RJLcBjwLHgOu6K3dg8OLvLcBZDK7a8codSZqgDHJ6dvV6ver3+9OehiSdVpLsr6rewrp/kStJDTH0Jakhhr4kNcTQl6SGGPqS1BBDX5IaYuhLUkMMfUlqiKEvSQ0x9CWpIYa+JDXE0Jekhhj6ktQQQ1+SGmLoS1JDDH1JaoihL0kNMfQlqSGGviQ1xNCXpIYY+pLUEENfkhpi6EtSQwx9SWqIoS9JDTH0Jakh66c9AUmadX94zbs4dvTvgZrCd38xN3xo36o9mit9SVrEIPD/jukEPsB/8K43XLNqj2boS9Iijn2lP+0pAP+xao9k6EvSYurItGewqsYO/STrknw6yYcX1N+WpJKcM1Tbk2Q+yeNJLh+qX5zkoW7fjUmyOm1I0hrJ2dOewapazkr/euCx4UKSLcAPAJ8fql0A7AIuBHYCNyVZ1+2+GdgNbO9uO0955pI0Aeuf25v2FIAXr9ojjRX6SeaA1wHvWbDrD4C38/9f4bgKuLWqjlbVE8A8sCPJecDGqrq3qgrYB1y90gYkaS1dv+8G1m+4DJjWiYnVvXpn3Es2380g3J/9PSfJlcAXqurBBWdpNgP3Dd0/0NW+2m0vrJ8gyW4GvxFw/vnnjzlFSVob1++7Abhh2tNYFUuu9JO8HjhUVfuHat8IvAN456ghI2q1SP3EYtXequpVVW/Tpk1LTVGSNKZxVvqXAlcmuQJ4HrAR+FNgG3B8lT8HfCrJDgYr+C1D4+eAg119bkRdkjQhS670q2pPVc1V1VYGL9B+vKp+tKrOraqtXf0A8Mqq+iJwJ7AryYYk2xi8YHt/VT0NHElySXfVzjXAHWvUlyRphFV/G4aqeiTJbcCjwDHguqp6ptt9LXALcBZwV3eTJE1IBhfSzK5er1f9/iz8RZwknT6S7K+qE643nfnQT3IY+Leh0jnAl6Y0nWlprWf7PbO11i9Mp+dvq6oTroSZ+dBfKEl/1LPXmay1nu33zNZavzBbPfveO5LUEENfkhpyOob+3mlPYApa69l+z2yt9Qsz1PNpd05fknTqTseVviTpFBn6ktSQmQn9JFuS/F2Sx5I8kuT6rv6hJA90tyeTPHCS8dcnebgb+5bJzn75Fun3oiT3df32u/czGjV+Z/chNfNJfnWys1++Vej3fUkOJXl4sjM/dSvp+WRjZ9kK+31ekvuTPNiN/fXJd7A8K/033R078sOp1lRVzcQNOI/B+/fA4C2c/xm4YMEx7wLeOWLsy4GHgW9k8NYSfwNsn3ZPp9IvcDfwg139CuDvR4xdB/wL8O3Ac4EHF/63mrXbSvrt9r0KeCXw8LR7mdD/4yV/HmbttsJ+A3xTt/0c4BPAJdPuaa36HXqMXwI+AHx4UvOemZV+VT1dVZ/qto8w+JSuZ99vv3uTtp8APjhi+MuA+6rqy1V1DPgH4IfXftanbpF+i8E7mQK8gNHvRLoDmK+qf62qrwC3Mvjwmpm1wn6pqn9kNT8degJW0vNSPw+zaIX9VlX9d3f3Od1tpq8yWem/6UU+nGptTfvZ8iTPflsZfATjxqHaq4D+SY5/GYNn2W9msNq/F/ijafdxKv12vXweeAr4AoM/pV54/I8B7xm6/9PAH0+7j7Xqd8G402alvxo9Lxw77T7Wsl8Gv8E+APw38NvT7mEC/f4FcDHwvbS40j8uyTcBfwm8par+a2jXTzJ6lU9VPQb8NnAP8FEGpzuOrfFUV8WIfq8F3lpVW4C3Au8dNWxEbaZXRcedYr+ntZX0vMjPw8w61X6r6pmquojBZ23sSPLySc15JU6l31EfTjUx036GXPDM9xzgY8AvLaivB/4dmBvzcX4D+Plp93Mq/QL/ydf/fiLAf40Y953Ax4bu7wH2TLuftep36NitnGYr/ZX0fLKfh1m+rfT/8dCYXwPeNu1+1qpf4DcZfA7Jk8AXgS8D75/EnGdmpd+ds38v8FhV/f6C3d8PfLaqDpw48tnx53Zfzwd+hJP8VjArFun3IPDqbvs1wOdGDP8ksD3JtiTPZfDhNneu5XxXaoX9npZW0vMSPw8zaYX9bkrywm77LLqf+bWd8cqspN8a/eFUb1zjKT/7zWfiBnw3g1MUn2FwXu8B4Ipu3y3Azy04/luBjwzd/ycGH9zyIPB90+7nVPvt6vu7Pj4BXHySfq9g8DrGvwDvmHY/E+j3g8DTwFcZrJDePO2e1rLnxX4eZvW2wn5fAXy6G/swI67Sm7XbSv9NDz3O9zLBc/q+DYMkNWRmTu9IktaeoS9JDTH0Jakhhr4kNcTQl6SGGPqS1BBDX5Ia8n8l+9ijlip/BwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "for gtx in ['gt1l','gt1r','gt2l','gt2r','gt3l','gt3r']:\n",
    "    flag = (D13[gtx]['h_ortho']>0)\n",
    "    y = plt.scatter(D13[gtx]['lat'][flag],D13[gtx]['h_ortho'])#[flag])\n",
    "plt.show()  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
